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two-dimensional  proof  of  Morawetz  that  such  flows  are  unstable 
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RESEARCH  INTO  TOPICS  IN  TRANSONIC  FLOW  THEORY 


(ONR  CONTRACT  N00014-83-C-0456 ) 


Summary 

The  contract  concerned  certain  research  topics  in  transonic 
flow,  mainly  addressing  the  questions  of  why  corrections  for 
entropy  and  viscosity  to  transonic  potential  theory  were  large  in 
certain  cases.  The  physical  implications  arising  from  the 
accurate  use  of  each  of  these  corrections  are  incompatible.  The 
investigation  concluded  that  the  main  cause  of  the  large 
corrections  is  due  to  the  fact  that  the  potential  equation  admits 
nonunique  solution  branches  and  that  the  corrections,  in  addition 
to  modeling  entropy  and  viscosity,  allow  the  potential  solution 
to  change  to  the  correct  branch.  In  view  of  this  hypothesis,  an 
analytic  study  into  the  causes  of  the  nonunique  potential 
solutions  was  initiated. 

The  study  into  the  causes  of  the  nonuniqueness  was  based  on 
an  integral  equation  formulation  and  the  following  results 
obtained . 

a.  The  nonuniqueness  only  occur  if  there  is  a  supersonic 
domain  in  the  flow. 

b.  A  nonunique  solution  does  not  have  the  same  supersonic 
domain  as  a  real  solution. 

c.  The  solution  often  can  be  made  unique  in  steady  flow  by 
using  a  mixed  conservative/nonconservative  algorithm; 
this  cannot  be  guaranteed  however. 


d.  There  is  no  way  to  differentiate  between  a  solution  that 
is  physically  correct  and  one  that  is  not  physically 
correct. 

An  additional  study  into  the  relative  magnitudes  of  the 
effect  of  vorticity  and  entropy  gradient  on  transonic  solutions 
indicated  that  both  are  equal  for  most  flow  regions.  However, 
the  shock  strength  is  controlled  only  by  the  entropy.  This  has 
relevance  to  the  "strong  shock"  potential  theories  in  which 
entropy  but  not  vorticity  is  retained  in  the  model. 

Finally,  an  investigation  into  the  stability  of  shock-free 
supercritical  flows  was  performed  for  three-dimensional  and  time- 
dependent  cases.  It  is  concluded  that  the  results  of  the  two- 
dimensional  proof  by  Morawetz  are  unstable  to  small  perturbations 
are  also  applicable  to  three-dimensional  and  unsteady  flows. 

A  detailed  report  on  the  work  is  attached  as  the  Appendix. 
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INTRODUCTION 


In  the  present  range  of  computational  fluid  dynamics  methods 
for  transonic  flow  there  is  a  hierarchy  of  equations  that  can  be 
solved.  For  practical  purposes  the  top  of  the  range  is  the 
Reynolds  averaged  Navier-Stokes  equations  with  some  form  of 
turbulence  model.  These  equations  will  model  transonic  flows 
with  viscous  effects,  strong  shock  waves  and  separated  flow;  the 
results  of  solving  these  equations  agree  quite  well  with 
experimental  data.  The  next  step  in  the  simplification  of  the 
flow  model  is  to  solve  the  Euler  equations  which  are  essentially 
the  Navier-Stokes  equations  with  the  viscous  effects 
eliminated.  The  elimination  of  the  diffusive  effects  leads  to  a 
loss  of  probable  uniqueness  in  the  problem  in  that  a  Kutta 
condition  at  the  trailing  edge  of  the  wing  may  be  necessary  to 
close  the  problem.  The  Euler  equations  will  model  flows  with 
strong  shock  waves  and  the  vorticity  caused  by  shocks.  The  final 
step  in  the  simplification  of  the  Navier-Stokes  equations  for 
transonic  flows  is  the  neglect  of  entropy  and  vorticity  in  the 
Euler  equations  which  leads  to  the  potential  equation,  or  for 
small  disturbances,  the  Transonic  Small  Disturbance  ( TSD ) 
equations.  These  potential  formulations  will  treat  any  attached 
flows  with  weak  shock  waves. 

In  certain  cases  the  potential  theory  can  be  coupled  with  a 
boundary  layer  model  to  treat  flows  with  mild  viscous  effects; 
for  example  see  Ref.  1.  This  modification  generally  weakens  the 
shock  wave  and  moves  it  upstream.  Another  recent  modification  of 
potential  theory  are  the  strong  shock  theories  pioneered  by  Nixon 
(Ref.  2)  and  Klopfer  and  Nixon  (Ref.  3),  which  use  a  potential 
theory  but  with  entropy  introduced  at  the  shock  wave.  The 
results  of  both  these  theories  are  considerably  better  than 
potential  theory  if  the  potential  shock  waves  are  strong.  In 
fact  the  correction  to  the  shock  locations  due  to  the  inclusion 
of  viscous  and  entropy  effects  is  nearly  the  same,  although  by 
different  mechanisms,  and  these  results  agree  relatively  well 


with  the  solutions  of  the  Navier-Stokes  equations.  The  strong 
shock  potential  theories  give  results  that  are  almost  identical 
with  the  Euler  Equations  solutions.  One  question  that  now  arises 
is  why  the  potential  theory  plus  boundary  layer  gives  such  good 
results  when  entropy  and  vorticity  effects  associated  with  the 
shock  wave  are  neglected;  this  implies  that  vorticity  and  entropy 
effects  are  negligible  relative  to  viscous  effects.  A  second 
question  is  why  does  the  strong  shock  theory  give  such  good 
results  if  compared  to  the  Euler  equations  results  since 
vorticity  is  neglected;  this  implies  that  vorticity  effects  are 
negligible  relative  to  entropy  effects  which  would  appear  to  be 
in  contradiction  to  Crocco's  Theorem.  Finally  a  comparison  of 
results  of  the  potential,  Euler  and  Navier-Stokes  equations 
indicate  that  viscous  effects  are  negligible  compared  with 
entropy  and  vorticity  effects.  The  above  noted  implications  are 
inconsistent  and  it  is  the  purpose  of  the  present  work  to  resolve 
at  least  some  of  these  inconsistencies. 

In  the  investigation  some  computational  results  led  to  the 
hypothesis  that  the  dominant  factor  in  the  large  entropy  and 
viscous  corrections  is  the  presence  of  the  nonunique  solutions 
first  noted  by  Steinhoff  and  Jameson  (Ref.  4)  and  further  studied 
by  Salas  et.  al.  (Ref.  5)  and  Nixon  (Ref.  6).  It  is  suggested 
therefore,  that  the  effects  of  viscosity,  entropy  and  vorticity 
are  relatively  small  compared  with  the  effects  of  a  nonunique 
solution.  Accordingly  a  study  into  the  causes  of  the 
nonuniqueness  was  initiated  and  the  following  conclusions 
reached . 

First,  it  can  be  shown  that  the  problem  does  not  occur  for 
subcritical  flows,  at  least  one  shock  wave  is  required  in  the 
flow.  Secondly,  the  solution  can  be  composed  of  a  real  lift 
component  and  a  ficticious  lift  component  due  to  the  nonunique 
behavior;  these  two  components  appear  as  a  pair  and  no  way  to 
decouple  them  could  be  found.  Thirdly  the  appearance  of  the 


nonunique  component  can  be  countered  in  steady  flow  by  the  use  of 
a  mixed  nonconservative/conservative  algorithm.  This  does  not 
guarantee  a  unique  solution  for  all  cases  however. 

A  second  part  of  the  study  concerned  the  relative  effects  of 
entropy  gradient  and  vorticity.  It  is  found  that  the  entropy 
gradient  and  vorticity  effects  are  the  same  magnitude  behind  the 
shock  but  that  there  is  an  additional,  perhaps  dominant,  effect 
of  entropy  in  determining  the  shock  strength  and  hence,  position. 

A  separate  study  concerns  the  nonexistence  proofs  for  shock 
free  flow  first  outlined  by  Morawetz  (Ref.  7)  and  alternately  by 
Nixon  and  Kerlick  (Ref.  8).  The  present  work  extends  the  theory 
of  Ref.  8  to  three  dimensional  steady  and  two  dimensional 
unsteady  flows. 

2.  VISCOUS  EFFECTS 

In  order  to  get  information  for  a  study  of  the  viscous 

effects  the  computer  code  GRUMFOIL  (Ref.  1)  was  used  to  compute 

the  flow  around  a  NACA  0012  airfoil  at  M  =  0.8,  a  =  0.75°.  It 

00 

proved  difficult  to  obtain  convergence  but  at  least  a  partially 

converged  result  (residual  of  ~  10-3)  was  obtained.  The  grid 

sequence  of  iterations  600/400/500  was  used;  the  fine  grid  is  160 

X  32.  The  results  are  shown  in  Figure  1.  The  same  problem  was 

then  solved  with  the  boundary  layer  omitted  and  this  result  is 

also  shown  in  Figure  1.  In  view  of  the  considerable  difference  a 

third  calculation  was  performed  namely  only  ten  iterations  on  the 

coarse  and  medium  grids  were  used  as  a  starting  solution  for  800 

iterations  on  the  fine  grid.  The  result  is  =ilso  shown  in  Figure 
*1 

±  • 

Because  of  the  disparity  in  results  the  same  case  was 
calculated  both  viscously  and  inviscidly  using  the  code  TAIR 
(Refs.  9,  10).  The  inviscid  results  are  compared  with  those  of 
GRUMFOIL  in  Figure  2,  the  viscous  calculations  in  Figure  3.  The 


TAIR  results  were  converged 


Apart  from  the  above  calculations  a  series  of  inviscid 
calculations  for  H#  =  0.8  and  different  angles  of  attack  were 
performed  using  GRUMFOIL  and  the  results  are  shown  in  Figure  4. 
One  set  of  results  is  for  the  usual  grid  sequencing  algorithm 
while  the  other  uses  10  iterations  on  the  coarse  and  medium  grids 
and  800  iterations  on  the  fine  grid. 

Although  the  fact  that  most  of  the  above  results  are  not 
fully  converged  raises  a  warning  flag  the  results  are  somewhat 
disquieting  to  the  author  since  it  appears  that  no  great  reliance 
can  be  placed  on  the  computations.  It  is  probable  that  the 
consistent  multiplicity  of  results  is  due  to  the  nonuniqueness  of 
the  potential  formulation  (Refs.  4-6)  and  that  the  viscous 
effects  are  in  fact  small,  of  the  order  of  that  shown  in  Figure 
3.  Consequently,  the  next  part  of  the  study  was  an  investigation 
of  the  uniqueness  of  the  potential  formulation. 

3.  NONUN I OUE  SOLUTIONS 

In  recent  years  multiple  solutions  to  the  numerical 
approximation  to  the  full  potential  equations  have  appeared  in 
the  literature  (Refs.  4  and  5).  Initially  the  phenomena  appeared 
in  computations  of  the  flow  over  a  symmetric  airfoil  at  zero 
angles  of  attack  where  two  lifting  solutions  were  present  in 
addition  to  the  expected  nonlifting  solution.  In  reference  4 
some  results,  for  a  nonsymmetric  airfoil,  a  RAE  2822  section,  are 
also  presented.  Steinhoff  and  Jameson  (Ref.  4)  suggested  that 
the  change  from  one  of  the  solutions  to  another  is  discontinuous 
and  noted  a  hysteresis  effect  indicating  that  the  lift 
coefficient  (C^)  depended  on  whether  the  angle  of  attack  (  a)  was 
increasing  or  decreasing.  More  recent  work  is  by  Salas  (Ref.  5) 
who  has  extended  the  computations  of  the  flows  considered  by 
Steinhoff  and  Jameson  (Refs.  4)  to  show  that  it  is  possible  to 
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construct  a  smooth  CL  -  a  curve  connecting  the  three  solutions 
for  a  symmetric  airfoil. 

The  investigations  noted  above  are  meticulously  performed 
and  are  essentially  numerical  experiments.  There  is  a  limited 
amount  of  understanding  that  can  be  gained  from  such  experiments 
and  consequently  a  more  analytic  technique  may  yield  additional 
information.  Futhermore,  although  the  numerical  results  are 
invaluable  they  do  not  exclude  the  possibility  that  the  multiple 
solutions  are  due  to  the  numerical  approximation  to  the 
differential  equation.  The  present  investigation  is  based  on  the 
integral  equation  formulation  (Ref.  11)  which  allows  some  degree 
of  insight  into  the  problem. 

The  transonic  integral  equation  method  of  reference  11  is 
only  applicable  to  the  transonic  small  disturbance  (TSD)  equation 
rather  than  the  full  potential  equation  (FPE)  that  is  used  in  the 
earlier  work.  Consequently  the  first  step  is  to  reproduce 
multiple  solutions  using  the  TSD  equation.  Once  these  solutions 
are  obtained  they  can  be  analyzed  using  the  ideas  of  the 
transonic  integral  equation  theory. 

3.1  Multiple  Solutions  for  Small  Disturbance  Theory 

Since  it  is  easiest  at  present  to  use  the  integral  equation 
theory  to  analyze  small  disturbance  theory  it  is  necessary  to 
reproduce  the  multiple  solutions  using  the  TSD  equation.  This  is 
achieved  by  using  the  computer  code  TSFOIL  (Ref.  12)  which  solves 
the  TSD  equation  using  the  conservative  Murman-Cole  algorithm  and 
grid  sequencing.  The  multiple  solutions  are  found  for  a 
symmetric  airfoil  at  zero  angle  of  attack  by  imposing  a  1°  angle 
of  attack  on  the  coarse  grid  solution  and  then  putting  the  angle 
of  attack  equal  to  zero  in  the  medium  and  fine  grid  operations. 

It  is  found  that  such  a  device  leads  to  multiple  solutions  over  a 
small  range  of  Mach  numbers.  Such  solutions  have  been  found  for 


a  11.8%  Joukowski  airfoil  and  a  NACA  0012  airfoil.  As  a  test  for 
convergence,  the  solution  for  the  Joukowski  airfoil  at 

*  0.85  was  converged  to  a  residual  of  10-7.  Krupp  scaling 
was  used  in  these  results;  the  default  Krupp  grid  is  used,  which 
has  grid  dimensions  of  77  x  56.  An  example  of  a  multiple 
solution  is  given  in  figure  5. 


3.2  Transonic  Integral  Equation  Theory 


Since  the  analysis  of  the  multiple  solutions  is  based  on  the 
integral  equation  theory  is  is  helpful  to  outline  the 
formulation . 

For  Krupp  scaling  the  TSD  equation  can  be  written  as 

(1  —  d>  )  <j>  +  6  =  0  (1) 

yx  yxx  vyy 


where  if  41  is  the  perturbation  velocity  potential  in  a  cartesian 
coordinate  system  (x,y)  then 


4>  =  k/B2  4 
x  =  x_ 
y  =  By 


where 


(2) 


k  =  (  y  +  1 )  M( 

2  2 
6  =  1  -  M 


1.75 


(3) 


If  u  =  4  ,  v  =  4  then 
'  x  ry 

(u,v)  are  given  by 


the  physical  perturbation  velocity  components 


u  =  (62/k)  u 
v  =  ( 8 3 /k )  v 


(4) 
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In  the  formulation  of  Equation  (1)  the  sonic  point  is  given  by 


u  -  ♦*-  1 


(5) 


The  boundary  conditions  for  Equation  (1)  are  that  4>  #  4>  +  0  on 

x  y 

the  far  field,  that  the  tangency  condition 


$  (x,±0)  =  y  (x,±0)  -  A 
y  s 


(6) 


is  satisfied  and  that  the  Kutta  condition  of  zero  velocity  jump  at  the 
trailing  edge  and  on  the  wake  is  satisfied.  In  Equation  (6) 


Yg(x,±0)  =  k/63Ys  (x , ±0 ) 


(7) 


where  y  =  Ys(x,±0)  denotes  the  geometry  of  the  airfoil  on  its  upper 
and  lower  surfaces,  respectively.  A  is  given  by 


A  =  k/e  a 


(8) 


where  a  is  the  angle  of  attack,  although  it  should  be  noted  that  in 
the  formulation  used  in  TSFOIL 

A  =  k/e3M  ~l/4a  (9) 

00 

The  basic  idea  of  the  integral  equation  method  is  to  use  Green's 
theorem  to  write  the  differential  equation.  Equation  (1),  and  its 
associated  boundary  conditions  in  integral  form.  A  detailed 
description  of  the  method  is  given  in  reference  13. 

For  y  f  0  the  integral  equation  is  given  by 

u  -  u2/2  =  u LS  +  uLA  +  IT(u)  (10) 


where 


ULS  (X'>° 


2tt  j0 


In 


(YsU,+0) 


YsU/-0)]  Kx(x,£;y,0)d£ 


(11) 


and  is  the  solution  of  Equation  (1)  without  the  nonlinear  terms, 


uLA(x,y)  /01Au  (c)Kn(xf C;y,0)d? 


(12) 


where 


Au(?)  =(u(C,+0)  -  u(C, -0)]/2  (13) 

The  kernel  function  K(x,£;y,ii)  is  given  by 

K(x,  £;y,n)  =  In  ( (x  -  U2  +  (y  -  n ) 2 ]  (14) 

The  integral  IT  is  given  by 

tT(j)  =  ~  17  W  Kcx(x,c;y/n)u2(e,n)ds  (15) 

The  domain  S  is  shown  in  figure  6. 


If  y  =  ±  0,  Equation  (10)  gives  (see  ref.  13)  only  the  symmetric 
part  of  the  solution  and  the  antisymmetric  part  is  given  by 


-A  +  i[Yg (x,+0)  +  Ys 


'  nn_  1  rl  [Au(0  -  Au  U)/2] 

«(x'  0,J"  7  Jo  - 


dC  +  I  (x) 
(16) 


where 


Au2 ( C )  =[u2(5,+0)  -  u2(C,-0)]/2 


(17) 


and 


I 

c 


(x) 


(x,c?0,n) 


[u2U,n)  -  u2(£)]ds 


where 


•  9 

u^U) 


u2( Sf+0) ,  n  >  0 
u2U,-0) ,  n  <  0 


(18) 


(19) 


This  version  of  the  theory  is  for  normal  shock  waves  and  the 
integrals  are  continuous. 

In  the  solution  of  the  integral  equations  the  circulation  is 
given  by  Equation  (16)  although  it  should  be  emphasized  that 
Equations  (10)  and  (16)  are  not  independent.  The  solution  of 
Equation  (16)  is  found  by  inverting  the  integral  equation  to  give 


su(xl  .  =  uLA  (.,01  -  i  (^)1/2  ll  ^ 


5  _  <1/2 


•  (rh) 


where  uLA  is  the  antisymmetric  solution  of  Equation  (1)  without 
the  nonlinear  terms  and  is  given  by 

,  ,  -a  +  i  (*'(;, +0)  +  y’(c,-oi] 

*  V  t^)172  /i  - 2"  S  „  S -  (t-^-f)17 


(x  -  5) 


irH)- 


The  inversion  procedure  invokes  the  Kutta  condition. 


3.3  Application  of  the  Integral  Equations  to  Multiple 
Solutions 


Consider  for  the  moment  the  case  of  a  symmetric  airfoil  at 
zero  angle  of  attack  in  this  case 


A  =  0  (22 

Y  ( X , +0  )  =  -Y  ( x , -0 ) 
s  s 

and  the  left  hand  side  of  Equation  (16)  is  zero.  A  manipulation 
of  the  integrals  in  Equation  (16)  leads  to  the  equation 

JSJ  K?y[u2U,n)  -  u2U,n)  -  2&uU)]cHdn  =  0  (23 

where  covers  half  of  the  domain  S.  Equation  (23)  can  be 
written  in  the  more  compact  form. 


JcJ  K_  f (5,n)d?dn  =  0 


where 


fU,ti)  =  fu2  (  £ ,  n)  -  u2U,-n)  -  2Au ( £ ) ] 


One  solution  of  this  integral  equation  is 

f  (  5  t  n)  -  0  (25; 

which  is  the  symmetric  solution  for  the  airfoil  problem.  A  non- 
symmetric  solution  can  be  obtained  if  there  is  one  or  more 
functions  f(£,n)  ?  0  that  satisfy  Equation  (24).  If  there  are 
such  solutions  then  a  multiple  solution  in  transonic  flow  can 
exist  if  the  equation 

\  [u2U,n)  -  u2U,-n)  -  2au(£)]  -  f(£,n)  (26) 

has  a  real  solution  for  u(£,ti). 


One  nonzero  form  of  f(£,n)  is 

1  /o  ?  0  <  £  <  1 

f(£,n)  -  U(1  '  5)1  C 

0  ;  £  <  o,  £  >  1 

where  C  is  an  arbitrary  constant.  This  is  the  form  of 
f(£,n)  that  is  equivalent  to  the  eigensolution  that  arises  in 
the  linear  airfoil  problem. 

Consider  new  Equations  (26)  and  (27)  on  n  =  0;  thus 


Au(£)  = 


(5(1  -  C )  1  x/2  (u  (  C  r  +0 )  +  u(£,-0)  -  2] 


It  can  be  seen  that  Au(£)  »  at  the  trailing  edge, 

£  “  1,  unless  the  average  velocity 

(u(£,+0)  +  u(£,-0)]/2 


-Wvi'V-:' >vvS>^ 


-10- 


behaves  like 


1  +  (1 


5) 


-P/2 


as  £  ♦  1,  where  p  >  1.  Hence,  unless  the  constant  C  is  zero, 
either  Au(£)  or  the  average  velocity  [u(£,+0)  +  u(£,-0)]/2  is 
infinite  at  the  trailing  edge.  The  condition  that  the  velocity 
is  finite  at  the  trailing  edge  is  thus  sufficient  to  eliminate 
the  eigensolution  by  putting  C  =  0.  This  is  a  slight  variation 
on  the  classic  use  of  the  Kutta  condition. 

This  integral  equation  is  valid  for  shock  waves  normal  to 
the  freestream;  if  the  shocks  are  not  normal  to  the  freestream  a 
modified  integral  equation  is  used  (ref.  14).  Since,  the 
formulation  changes  are  negligible  the  above  set  of  equations 
will  be  used  in  the  subsequent  algorithm  for  clarity. 

3.4  Number  and  Nature  of  Eigensolut ions 


It  is  desirable  to  estimate  the  number  and  the  nature  of  the 
eigensolutions ,  f(£,n)  of  equation  (24).  Consider  equations 

(16)  and  (18).  Using  the  inversion  procedure  that  leads  to 
equation  (20),  equations  (16),  (18)  can  be  written  in  the  form 


K?y  [(u2U,n)—  u2U,-n) 


-2 Au ( £ )  -21  (£)]d£dn  =  0 
C  L 


(29) 


where 


JL 

IT 


(1  -  x } 1/2 


a  -A  +  i  <*»<«>  +  V«>) 

O  (X  -  £ ) 


X/2d£ 

(30) 


and  u(£,n)  is  a  transonic  solution  subject  to  the  arbitrary 
boundary  conditions  in  equation  (6).  Since  the  boundary 
condition  is  arbitrary  it  follows  that  an  infinite  number  of 
solutions  u(£,n)  exist. 


Equation  (29)  is  identical  in  form  to  equation  (24)  with 
f ( ?, n)  given  by. 

f(?,n)  *  u2U,n)  -  u2 ( 5 ,-n )  -  2au(?)  -21  (?)  *  0  (31 

CL 

Since  there  are  an  infinite  number  of  transonic  solutions  it 
follows  that  there  are  an  infinite  number  of  eigensolutions 
f(?,n)  .  A  justification  of  a  nonzero  f(?,n)  is  given  in  the 
Appendix. 

3.5  Symmetric  and  Antisymmetric  Integral  Equations 
Let  u  (?,n)  and  Au(?,n)  be  defined  by 


u  ( ? f  n )  =  (u(?,n)  +  u(?,-n)] 


(32a 


Au( ?,n)=  ~  (u(?,n)  -  u  (?,-n)} 


(32b 


Equation  (10)  can  be  manipulated  to  give  the  following  symmetric 
and  asymmetric  parts. 


u(x,y)  -  i  [ u 2 ( x , y )  +  Au2(x,y)]  =  uL~  ~  fg  J K  [u2(?,n) 


+Au  (?, n) ]  d?dn 

Au ( x , y )  -  Au(x,y)u(x,y)  =  -  /g  /K5x[u2(?,n)  -u2(?,n) 

-2au ( ? ,o) ) d?dn 

These  equations  will  be  used  in  the  following  sections 


3.6  First  Compatibility  Condition 

If  in  equation  (24),  f(?,n)  is  assumed  known  then  equation 

(23)  gives 


-  [u‘(?,n)  -  u  ( £,n)-2Au(  5,0)]  =  Au(e,n)uU,n)  -  Au(e,0)  =fU,n) 
2 


If  u(£,ri)  is  known  then  this  is  an  equation  for  Au(£,n)  in  the 
flow  field  and  gives 


AUq  - 


uo-1 


where 


fQ  =  f (5 ,0)  etc. ; 


also 


f  +  Au 


Au  = 


2  =  i  [£  -  -£2_] 


In  the  normal  transonic  solution  equation  (16)  only  gives  the 
value  of  Au ( C , 0 )  or  auq,  the  value  of  Au(x,y)  being  found  from 
equation  (34).  Hence  for  a  nonunique  solution  to  exist  Au,  as 
defined  by  equations  (36)  and  (37)  must  be  compatible  with 
equation  ( 34 ) . 

Substitution  of  equations  (35),  into  equation  (34)  gives 

AU-Auq  =  f  ~  f  d^dTl  (38) 

or,  using  equations  (36),  (37) 

-  (f  -  fo  f  -  i-  /  / K_  f  d?dn  (39) 

u  ( l-uQ )  411  bl 

An  alternative  form  of  equation  (39)  is  the  differential  equation 


where 


F(x,y ) 


J* 

/X 


1 

u  (5»y) 

[  AuU»y) 


f U,y)  -  f  (5#o) 
-  AuU,0)}  d£ 


[l-u(s,y )) 
[l-u( £,0)] 


dC 


(42 

(43 


The  far  field  boundary  condition  is  that  Fx,Fy  0. 

The  basic  differential  equation,  (1)  can  be  decoupled  into 
symmetric  and  asymmetric  parts.  The  asymmetric  part  is 


A  d>  + 

Txx 


A  4>  =  (  AUU  ) 

yy  x 


(44 


where 


A<J>(x,y)  =  j  U(x,y)  -  4>(x,-y)] 


(45) 


the  boundary  condition  is 


A  <|>y  ( x  , o )  ■  -A  +  1  [  Yg  (  x  ,  +0  )  +  Yg  (x  ,-0  )  ] 


(46) 


The  symmetric  part  is 


—  1  —9  9 

(,  +  <k  =4  (u  +  Au  1 

*xx  *yy  2  x 


( 47  • 


with 


^  ( X , o )  =  \  [Yg(x,+0)  -  Yg  (x,-0 ) ] 


(48 


Now  equation  (40)  can  be  modified  to  give 


where 


G<x,y)  *  /  x  -  g(?,y)  <3?  (50) 

u  U,y) 

and 

g(5,y)  »  f(C,y)  — —  (51) 

l-u  (5,0) 

with  the  boundary  condition 

Gy(x,±0)  =  0  (52) 

If  G(5,y)  is  identified  with  a  lifting  term  A <f>  then 
equations  (47)  and  (49)  can  be  added  to  give 

U  +  G)xx  +  {*  +  G)yy  =  \  C(u  +  5)2]x  (53) 

u 

1  ( x,  ±0  )  +  G  (x,±0)  C  Yg  ( X , +0 )  -Yg  ( x ,  -0  )  ]  (54) 

where  Au  on  the  right  hand  side  of  equation  (47)  has  been 
identified  as  g/  u  .  Thus  a  ficticious  value,  G,  can  be  added  to 
a  purely  symmetric  problem,  denoted  by  to  give  lift;  the 

boundary  condition  is  not  affected.  This  is  the  mechanism  of  the 
appearance  of  the  nonunique  solutions.  Equations  (53)  (54)  are 
similar  in  form  to  equations  (1)  (6)  if  4>  ±s  identified  with 
"$  +  G.  The  numerical  algorithm  will  solve  equations  (53)  in  an 
identical  manner  to  the  real  solution. 

If  there  is  real  lift,  that  is  a  value  that  is  calculated 
with  no  eigensolutions,  then  the  lifting  solution  of  equations 
(44)  can  be  added  to  equation  (53)  to  give 

("♦  +  A  4>-*-G )  +  ("$  +  A<|>+G)  =  i  [  ( u  +  Au+C[)2]  (55) 


where  denotes  a  real  lifting  component.  In  this  case  au  on 
the  right  hand  side  of  equation  (47)  is  identified  with 
Cau  +  g/-3. 

The  boundary  condition  is 

1  (x,±0)  +  A**(x,±0)  +  G  ( X , ±0 )  -  Y ' ( x , ±0 )  -A  (56) 

y  y  y  s 

Again  a  ficticious  lifting  component  G(x,y)  is  added  to  the 
equation  without  a  change  in  the  boundary  conditions.  If  G  is 
present  the  numerical  algorithm  will  solve  equations  (55)  in  an 
identical  manner  as  the  real  solution.  It  is  not  possible  to 
decouple  the  real  and  ficticious  components. 

3.7  Second  Compatibility  Condition 

In  the  preceeding  analysis  it  has  been  assumed  known.  For  a 
solution  to  exist  it  also  must  be  a  solution  of  equation  (33). 

Equation  (39)  can  be  written  as 


where  g  is  given  by  equation  (51).  It  follows  from  equation  (57) 
that 


u 


g  +  if 


it  also  follows  from  equation  (58)  that 


A u  =  g/-  =  g  +  If 


(58) 


(59) 


From  equation  (36) 


■.»  v  V  ^  v 


i„ 

-r4 


I 


!_ 


•  .*> 
•’A 


L- 


and  if  the  limit  as  y  ♦  0  is  taken  of  equation  (33)  and  if 
equations  (58),  (59)  are  used  then 


-2 

u 


uo  -  r  -  “ls„  +  5  <37- r>2  -  k  *is  Is,  /Kcx  [Ti«7T2  +  <«+If )2]d«d" 

f  (60) 

where  uLg  is  the  value  of  uLg  at  y  =  0. 

Since  f  is  assumed  known  and  uLSq  is  known,  equation  (60)  is  an 
equation  for  vJo  .  Hence  g  can  be  found. 

Using  equations  (58),  (59)  to  eliminate  u  and  au  in 
equation  (33)  it  follows  that  a  nonunique  solution  only  exists  if 
there  is  a  solution  f,  to  the  equation 

g+If  ~  2  (g+If  J  =  ULS  +  2  (9+If)  "  4^  ^KxC*-  (g+i  _)'2  +  (^+If)  3<Hdn 


(61) 


This  solution  must  also  satisfy  equation  (24). 


If  there  is  real  lift  due  to  the  boundary  condition  of  equation 
(6)  then  equation  (24)  gets  replaced  by 


/fi.  /  KCy  Cf+Au*u*-  au*  -ln  ]d£dn 


o  cT 


(62) 


*  _  * 

where  Ic  is  given  by  equation  (30)  and  Au  ,u  denote  a  lifting 

L  i, 

and  a  symmetric  component  of  some  solution  u  that  satisfies  the 


boundary  condition  represented  by  Ic  , 

*  L 

u  cannot  occur  at  the  same  location  as  in  the  real  solution,  u 


Since  the  shock  waves  in 

* 


cannot  be  the  real  solution. 


If  f  is  defined  by 


*_*  * 
Au  u  -  Au 

o 


(63) 


i". 


-1 


"•-'H 


-17- 


then  it  follows  that  equation  (26)  is  replaced  by 


Auu  -  auq  »  f  +  f*  (64) 

The  analysis  in  section  3.6  and  the  present  section  is  unchanged 
if  f  is  replaced  by  f  +  f*  since  f  always  occurs  as  part  of  the 
term  f  +  f* ;  in  this  case  f  +  f*  must  satisfy  equation  (62) 
rather  than  equation  (24). 

3.8  General  Remarks  on  Integral  Equations 


Consider,  for  example,  equation  (10),  which  can  be  written  as 


u  -  u2/2  -  +  IT  (u2) 


(65) 


The  integral  IT  can  be  discretized  to  give 


u 


i 


^2 


N 

+  I 

j=l 


u. 


(66) 


where  A^j  is  an  influence  coefficient  and  N  is  the  number  of 
discrete  elements  used  in  the  evaluation  of  the  integral. 


The  algebraic  set  of  equation,  (66)  is  quadratic  and  has,  in 
general,  2N  solutions.  As  N  -*■  »  there  are  an  infinite  number  of 
solutions.  However,  if,  during  the  solution  of  u^  say,  only  one 
root  is  taken  then  the  number  of  solutions  is  1N  or  in  other 
words  a  unique  solution  is  obtained.  Thus  equation  (66)  will 
give  a  unique  solution  if  one  of  the  two  choices  of  u^  is 
eliminated.  A  multiplicity  of  solutions  occurs  if  a  choice  is 
allowed. 


In  the  general  transonic  solution  a  change  of  sign  occurs  at 
the  sonic  line  and  the  shock  wave  (Reference  13).  The  correct 
solution  is  that  choice  of  root,  or  sign,  that  eliminates 
expansion  shocks.  A  change  in  the  location  of  the  sign  change 
indicates  an  alteration  in  the  size  of  the  supersonic  domain. 


Equation  (61)  can  be  written  as 


'9+1, 


+  V  9)  -  j  (g+f 


+  if+  g)  =  uLS+  -  4it  fs  f  Kx  Cj-+~ 


+  g+  Ij]  d£dn 


(67) 


where 


LA 


-fo/(1-uo) 
(x-5)2  + 


<3? 


In  equation  (67)  there  is  the  possibility  of  two  solutions  for 
[j +  If+  g];  this  is  equivalent  to  u.  If  only  the  real 
solution  is  considered,  that  is,  u*,  the  value  of  u  corresponding 
to  f*,  then  equation  (61)  will  have  a  specified  combination  of 

sign  changes.  If  the  location  of  these  sign  changes  does  not 

★ 

alter,  then,  since  u  must  be  the  only  solution  of  the  system,  no 
ficticious  term  can  arise.  A  fictitious  term  u  can  arise  if 
there  is  an  alteration  of  the  sign  changes,  that  is,  the  size  of 
the  supersonic  domain  changes.  Thus  in  general,  a  nonunique 
olution  will  change  the  location  of  the  shock  wave. 


In  a  subsonic  flow  only  one  root,  the  subsonic  (negative) 
root  is  chosen  and,  since  no  sign  change  is  possible,  it  follows 
that  a  nonunique  solution  cannot  occur.  Thus  a  nonunique 
solution  can  only  occur  if  there  is  at  least  one  shcck  wave  in 
the  flow.  It  should  be  noted  that  there  is  no  information  about 
the  number  of  nonunique  solutions  since  this  depends  on  the 
number  of  possible  shock  locations. 

3.9  Non-Conservative  Algorithms 


A  non-conservative  algorithm,  such  as  that  of  Murman  and 
Cole.  (ref.  15)  adds  source  terms  at  the  shock  waves.  Hence  the 
algorithm  is  solving  conservatively  a  differential  equation  of 
the  form. 


where  o. (y )  is  a  source  term  at  the  shock  location  x_  and  H  (  ) 
1  si 

is  the  step  function.  The  strength  of  the  source  term  is  not 

known  explicitly  in  nonconservative  algorithms. 

If  the  asymmetric  part  of  equation  (68)  is  taken  then 

A4>xx+  Myy  ■  Uuu)x+  A[0iH(x-xs  )]x  (69) 

If  equation  (69)  is  added  to  equation  (49)  the  following  result 
is  obtained. 

U$+G)xx  +  (A^+G)^  =  (Auu)x  +  (g  +  A[oiH(x-xs  ])x  (70) 

If  a  ficticious  positive  component  of  lift  starts  to  appear 
in  a  solution  g  will  change  from  zero  to  a  predominantly  positive 
quantity,  since  g  is  equivalent  to  (Auu).  This  introduction  of 
positive  lift  will  change  the  location  of  the  shock  waves  such 
that  an  upper  surface  shock  moves  aft  while  a  lower  surface  shock 
moves  forward.  It  follows  that  the  source  term  introduced  by  the 
shock  on  the  lower  surface  now  acts  at  a  further  forward  location 
while  the  countering  source  on  the  upper  surface  moves  aft. 
Consequently,  there  is  a  considerable  region  between  the  shocks 
for  which 

A[ o. H (x-x  )]  <  0 
1  si 

This  counters  the  positive  g  and  reduces  the  tendency  of  the 
nonunique  solution  to  appear.  However,  it  is  necessary  to  point 
out  that  the  nonconservative  source  error  must  be  large  enough  to 
counter  the  g  terms.  The  nonconservative  algorithm  only  inhibits 
a  nonunique  solution  from  appearing;  it  will  not  necessarily 
remove  an  existing  error.  A  similar  analysis  can  be  performed  if 
the  fictitious  lift  is  negative,.  This  hypothesis  was  tested  by 


using  the  computer  code  TSFOIL  (Ref.  12)  which  has  both 
conservative  and  nonconservative  algorithms.  A  composite 
algorithm,  given  by 

Algorithm  =  A  (conservative)  +  (1  -  A)  (nonconservative)  (71) 

was  used.  The  parameter  ,  A,  was  taken  to  be  given  by 

A  =  -e|[CL(t+At)  -  CL(t)]|  +  1  (72) 

Where  t  is  the  artificial  iteration  time  and  At  is  the  iteration 
step.  At  convergence  A  =  1  and  the  solution  is  conservative. 

It  was  found  that  for  the  nonunique  solutions  to  be  avoided  in  a 
computation  of  the  example  in  Figure  (5)  e[C  (t+At)  -  C  (t)]  had 
to  be  of  order  unity  during  most  of  iteration.  Smaller  values 
did  not  inhibit  nonuniqueness.  Hence  it  is  suggested  that  the 
"classic"  degree  of  nonconservative  algorithm  is  necessary  to 
stop  the  appearance  of  nonunique  solutions. 

It  is  of  interest  to  note  that  if  a  stabilizing  term  of  the 
type  used  in  TSFOIL  namely  e$xt  is  added  to  equation  (1),  where 
t  is  an  artificial  iteration  time  and  e  is  a  parameter,  then  an 
analysis  similar  to  that  given  in  this  section  shows  that  this 
term  will  assist  the  formation  of  a  nonunique  solution.  This  is 
possibly  the  reason  why  the  real  solution  cannot  sometimes  be 
computed . 

3.10  Summary 

The  analysis  given  in  this  section  can  be  summarized  as 
follows . 

The  TSD  equation  can  admit  eigensolutions  that  satisfy  all 
of  the  boundary  conditions  generally  found  in  such  problems; 
there  are  an  infinite  number  of  these  eigensolutions.  For  a 


nonunique  solution  to  exist  certain  compatibility  requirements 
must  be  met.  These  eigensolutions  provide  lift  and  act  like  an 
additional  asymmetric  source  term  in  the  flow  field.  If  real 
lift  is  present  the  fictitious  component  provided  by  the 
nonuniqueness  appears  as  a  simple  additive  term  to  the  real 
lifting  component;  there  is  no  way  to  distinguish  or  uncouple  the 
two  components  in  a  given  numerical  solution. 

The  nonuniqueness  does  not  appear  for  subcritical  flows  and 
if,  in  a  supercritical  flow,  it  does  appear  the  location  of  the 
shock  waves  changes  from  their  real  location.  The  nonuniqueness 
may  be  removed  in  some  steady  cases  by  using  a  nonconservative 
conservative  algorithm  which  inhibits  the  appearance  of  such 
solutions.  However,  there  is  no  guarantee  that  this  will  work  in 
all  cases.  This  partial  cure  only  works  for  steady  flows;  it  is 
possible  that  for  unsteady  flows  that  a  "strong  shock"  theory  of 
the  type  advocated  in  Ref.  3  could  remove  any  nonuniqueness. 

Since  the  behavior  of  the  nonunique  solutions  is  identical 
to  that  of  the  real  solution  the  question  arises  as  to  whether 
these  solutions  are  physically  realizable.  An  analysis  of  the 
Navier-Stokes  or,  possibly,  the  Euler  equations  is  necessary  to 
determine  this. 

4.0  ENTROPY  EFFECTS 

The  motivation  for  examining  the  magnitude  of  the  entropy 
effects  is  the  good  agreement  of  the  results  of  the  strong  shock 
theory  of  Klopfer  and  Nixon  (Ref.  3)  with  results  of  the  Euler 
equations.  This  good  agreement  implies  that  entropy  effects  are 
large  compared  to  vorticity  effects.  This  implication  is  studied 
in  this  section.  The  analysis  technique  will  be  based  on  an 
integral  equation  for  the  Euler  equations  developed  by  Nixon  and 
Liu  (Ref.  16). 
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4.1  An  Integral  Equation  Formulation  for  the  Euler  Equations 


For  the  present  purpose  only  two-dimensional  flow  is 
considered.  It  is  not  necessary  to  consider  the  full  set  of  the 
Euler  equations  at  present  and  here  only  the  conservation  of  mass 
equation 
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and  the  vorticity  equation 
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are  considered,  where  p  is  density  and  u  is  vorticity. 


Equation  (73)  can  be  satisfied  by  a  stream  function  where 
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Substitution  of  equation  (75)  into  equation  (74)  gives 
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where 
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(77) 


Thin  airfoil  boundary  conditions,  namely 


V  ( x ,  ±0 )  =  -[a  +  Y  ( x ,  ±0  )  ]U 
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(78) 


are  assumed  to  be  applicable,  and  it  is  also  assumed  that 
freestream  conditions  apply  ahead  of  and  above  and  below  the 
wing.  At  the  downstream  boundary  constant  (freestream)  pressure 
is  imposed. 


Following  the  analysis  of  Nixon  and  Liu  (Ref.  16)  equation 
(76)  can  be  put  in  integral  form?  thus 
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where 


K  =  i  In  {(x-c)2+  (y-n)2)  (80) 

and 

Ko  -  \  In  {(x-c)2  +  y2}  (81) 

Dif ferentation  of  equation  (79)  with  respect  to  y  gives 

*y  '  ‘  h  /o  U*„Kc=y  -  Koyna»)d{  +  Ti  IJ  Ky£ds  (82) 

By  integrating  the  second  term  in  the  line  integral  by  parts  with 
respect  to  £  and  the  field  integral  by  parts  with  respect  to 
n,  the  following  equation  is  obtained. 
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F(x,y)  =  /yf (x, n)dn 
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Dif ferentation  of  equation  (79)  with  respect  to  x  gives 
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Performing  the  same  integration  by  parts  as  above  and  letting 
y  -*■  0  gives 


where  the  subscript  "0"  denotes  a  value  at  y*=0 .  In  equation  (86) 
the  left  hand  side  is  determined  from  the  boundary  condition, 
equation  (78)  p  is  known. 


Consider  the  behavior  of  equation  (83)  as  x  ».  The 
boundary  condition  for  a  flow  variable  P  as  x  ♦  is 
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Now,  since  the  kernel  K  vanishes  as  |x  -  g|  ♦  •  it  follows 
that  as  x  -*■  »  only  the  value  of  F(  5,  n)  relatively  close  to  x 
contributes  to  the  integral.  As  -r—  ♦  0  as  x  +  *  is  a  boundary 
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condition  F(?,n)  can  be  replaced  by  F(x,n)»  Thus 
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This  integral  is  only  finite  if 


Lim  U(pU)  -  AF]  =  0 
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In  Reference  (16)  the  use  of  the  boundary  condition  of  freestream 
pressure  at  x  -*■  «  and  equation  (87)  indicates  that  for  symmetric 
flows 


F 

Hence  equation  (87)  is  consistent  with  the  result  of  Reference 
(16). 

A  final  point  is  that,  from  equation  (86), 
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If  the  orthogonal  coordinate  system  (x,y)  is  chosen  such 
that  the  x  axis  aligns  with  the  slit  representing  the  airfoil  and 
the  wake  it  follows  that  equation  (88)  is  satisfied,  since  there 
is  no  flow  through  the  wake  and  hence 

( pV )  =0  on  y  =  0  (92) 

o 

4.2  Small  Disturbance  Versions  of  The  Euler  Equations  in  the 
Wake 

It  is  convenient  now  to  make  a  small  disturbance  assumption 
in  the  wake.  Although  the  velocity  over  most  of  the  airfoil  may 
be  large  and  these  determine  the  entropy  and  vorticity 
distributions  it  is  plausible  to  assume  that  at  the  trailing  edge 
of  the  airfoil  and  on  the  wake  a  small  disturbance  form  for  the 
velocity  can  model  the  flow.  Hence,  in  the  region  of  the 
trailing  edge  and  the  wake  it  is  assumed  that 

p  =  pm(1  +  p  ) 

U  =  (1  +  u)  (93) 

V  =  v 

and  it  is  assumed  that 

|p  | ,  | u | ,  jvj,  <  <  1  (94) 

Furthermore,  it  is  assumed  that 

0  j^|  ~  0  (u)  in  the  wake  (95) 

Where  S  is  entropy  and  R  is  the  gas  constant.  Along  the  wake 
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From  the  energy  equation  (see  Ref.  16) 
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where  Cp  is  the  specific  heat  at  constant  pressure  and  am  is  the 
sonic  speed  at  freestream  conditions. 


From  equation  (97) 
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from  Crocco's  theorem 
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Where  cv  is  the  specific  heat  at  constant  volume 
Using  the  gas  relation 
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Equation  (98)  can  be  written  as 
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Using  equations  (99)  and  (101)  in  equation  (97)  gives 
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(100) 


(101) 


(102) 


If  the  small  disturbance  assumption  is  valid  in  the  wake  and 
linear  approximations  are  valid  equation  (102)  is  approximated  by 
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and  hence 
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since  there  is  zero  pressure  jump  across  the  wake.  In  this  case 
the  limits  of  integration  in  the  line  integral  of  equation  (86) 
are  replaced  by  0,  1  and  an  explicit  wake  calculation  is  not 
necessary  however,  the  wake  geometry  does  appear  through  the 
field  integral.  This  is  a  similar  situation  to  the  behavior  of  a 
solution  of  the  full  potential  equation  with  nonlinear  boundary 
conditions . 


If  the  small  disturbance  assumption  of  equations  (93)  -  (95) 
are  used,  then,  in  the  wake 
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where  the  pressure  relation 
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has  been  expanded  for  small  u  and  S  .  Using  equation  (104)  in 
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equation  (86)  yields 
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This  equation  can  be  inverted  to  give 
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where 
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In  the  inversion  use  has  been  made  of  the  fact  that  at  the 
trailing  edge 
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This  is  a  Kutta  condition  and  is  necessary  to  ensure  a 
finite  value  of  [A(pV)  -  AF]  at  the  trailing  edge.  However,  if 
the  solution  algorithm  for  the  Euler  equations  automatically 
captures  a  contact  discontinuity  at  the  trailing  edge  across 
which  Ap  is  zero  then  such  a  condition  is  superfluous.  The 
necessity  of  a  Kutta  condition  is  determined  by  whether  or  not 
the  algorithm  captures  a  contact  discontinuity  at  the  trailing 
edge . 

It  is  instructive  to  examine  the  wake  problem  if  a  "strong 
shock"  potential  theory  (e.g.  Ref.  3)  is  used. 


The  theory  of  Ref.  3  neglects  vorticity  while  retaining 
entropy.  If  the  small  disturbance  approximation  is  used  in  the 
wake,  equation  (103)  is  replaced  by 
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and 
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A  (  pU ) -  aF  =  -  y— ^  AS  on  the  wake  (112 


Substitution  of  equation  (112)  into  equation  (107)  gives 
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Since  the  entropy  jump  is  constant  along  the  wake  and  is 
determined  by  the  shock  strength,  the  second  term  on  the  right 
hand  side  acts  like  a  boundary  condition.  Since  AS  is  constant 


this  integral  is  infinite  and  is  thus  inappropriate.  In  the 
numerical  examples  in  Ref.  3  the  downstream  boundary  is  at  a 
finite  distance  from  the  airfoil.  In  such  cases  the  entropy 
integral  is  finite  but  is  dependent  on  the  location  of  the 
boundary.  It  appears  that  a  consistent  strong  shock  theory  is 
not  possible?  for  results  that  are  independent  of  the 
computational  boundary  AS  should  be  zero  across  the  wake. 


It  should  be  noted  that  if  AS  is  positive  then  the  entropy 
integral  acts  to  produce  a  negative  camber  and  thus  reduces  the 
lift.  This  accounts,  perhaps,  for  the  excellent  agreement  of  the 
method  of  Ref.  3  with  the  results  of  the  Euler  equations. 


Equation  (113)  is  similar  to  one  used  in  aeroelasticity  and 
can  be  inverted  to  give 
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where 
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and  A  is  given  by  equation  (109).  In  the  inversion  a  Kutta 
condition  ensuring  that 


CA(pU)  -  A F]x-1  =0  (116) 

must  be  enforced.  If  the  small  disturbance  formulation  is 
applicable  at  the  trailing  edge,  then,  neglecting  vorticity, 

A  (  pU )  -  AF  =  pb  Au 


(117) 


and  thus  the  Kutta  condition  enforces  a  zero  velocity  jump  at  the 
trailing  edge  for  a  finite  solution  rather  than  a  zefo  pressure 
jump  as  is  the  case  with  the  Euler  equations.  In  other  words, 
the  strong  shock  theory  requires  a  Kutta  condition  to  make  the 
velocity  jump  finite  at  the  trailing  edge  and  this  condition 
disallows  the  possibility  of  making  the  pressure  jump  zero  unless 
as  «  0. 


4.3  Shock  Location 

The  most  noticeable  difference  between  potential  and  Euler 
equation  solutions  is  the  shock  location  consequently  it  is  of 
interest  to  investigate  the  mechanism  that  fixes  the  shock 
location.  In  the  following  analysis  only  non-lifting  flows  are 
considered. 

The  integral  equation  for  non-lifting  flows  can  be  derived 
from  equation  (83)  and  gives 

pU  =  I-  £1( pV)Kox  a«  -  k  I  I  KnyFds  die) 
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The  shock  location  appears  only  through  the  field  integral  over  S 
so  that  the  location  of  the  shock  depends  on  the  magnitude  of 
F.  The  analysis  can  be  simplified  if  the  field  integral  is 
integrated  by  parts  to  give 

JSJ  K^yFdS  =  -  JgJ  KyfdS  (119) 

and  f  is  given  by  equation  (77) 


The  vorticity  on  the  airfoil  surface  is  given  by  Crocco's  theorem 
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The  density  is  given  by 
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If  the  function  f  is  expanded  about  the  freestreara  conditions 
then 

f  -  isr  §  ~  ^  §  -  p.  M!  u- 17  +  second  order  terms  (123> 

The  first  term  is  an  approximation  to  the  vorticity  term 
while  the  remaining  first  order  terms  are  an  approximation  to 
PyU  .  It  can  be  seen  that  the  vorticity  term  is  small  compared 
to  the  entropy  effect  in  the  density  only  if 
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Note  that  both  the  vorticity  and  entropy  contributions  have  the 
same  sign.  Hence  the  inclusion  of  only  one  component  is  probably 
an  improvement  on  the  potential  theory. 


This  condition  is  not  satisfied  generally  in  transonic  flows. 
Hence,  in  this  case  vorticity  effects  and  entropy  gradient 
effects  are  the  same  order  of  magnitude. 

There  is  another  effect  of  entropy  on  a  transonic  solution, 
namely  the  effect  on  the  shock  jump  condition.  This  weakens  the 
shock  wave  and  moves  the  location  forward.  The  magnitude  of  this 
effect  is  difficult  to  estimate  theoretically.  As  a  general 
statement,  vorticity  cannot  be  neglected  relative  to  entropy;  if 
this  is  done  only  a  part  of  the  necessary  overall  correction  is 
retained.  However,  the  entropy  correction  alone  is  probably 
better  than  no  correction. 


4 . 4  Summary 

In  this  section  it  has  been  shown  that  vorticity  cannot  be 
neglected  if  entropy  is  retained  in  the  flow  equations.  However, 
if  only  the  entropy  effect  is  retained  this  is  probably  an 
improvement  on  potential  theory.  It  has  also  been  shown  that  if 
the  entropy  jump  across  the  wake  is  retained  in  a  strong  shock 
theory  the  solution  will  be  dependent  on  the  size  of  the 
computational  boundary;  it  is  also  impossible  to  have  both  a 
finite  value  of  aU  and  a  zero  pressure  jump  at  the  trailing 
edge. 

5.0  PERTUBATIONS  OF  A  SHOCK  FREE  SUPERCRITICAL  FLOW 

In  the  early  days  of  transonic  flow  research  the  question  of 
the  existence  of  a  shock  free  supercritical  flow  was  posed.  The 
question  was  answered  in  a  series  of  papers  by  Morawetz  (7,17) 
who  used  the  hodograph  transformed  equations  to  show  that  a  shock 
free  supercritical  flow  is  unstable  to  arbitrary  small 
perturbations.  This  instability  was  regarded  as  proof  that  such 
flows  could  not  exist.  However,  several  years  later  the  work  of 
Nieuwland  (18),  also  using  the  hodograph  transformation, 
indicated  that  shock  free  supercritical  flows  could  exist.  These 
two  apparently  conflicting  results  were  reconciled  by  realizing 
that  shock  free  supercritical  flows  could  exist  but  were  unstable 
to  a  small  perturbation. 

The  hodograph  transformation  used  in  the  analysis  of 
Morawetz  cannot  be  used  for  either  three  dimensional  or  unsteady 
flows.  In  the  present  note  the  transonic  integral  equation 
method  (11,13)  is  used  to  examine  the  perturbation  of  a  shock 
free  supercritical  flow.  The  analysis  is  a  minor  extension  of 
that  given  in  Reference  (8). 


5 . 1  Some  Comments  on  Integral  Equations 


Since  the  analysis  is  to  be  based  on  the  transonic  integral 
equation  method  some  basic  points  are  outlined  below.  Full 
details  of  the  particular  formulation  can  be  found  in  Reference 
13. 


The  three-dimensional,  low  frequency  transonic  small 
disturbance  equation  is 
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z  is  normal  to  the  wing. 

This  can  be  written  in  integral  form  to  give  (Ref.  14,  19). 
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where 
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the  kernel  K  is  given  by 
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K(x, n;z, c)  =  l/C(x-£)2+  (y  -  n)2+  (z-;)2]2  (129) 

A  is  the  planform  of  th  wing  and  wake,  S  is  the  flow  domain  after 
removing  the  wing,  wake,  shock  waves  and  the  field  point  (x,y). 
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For  a  shock  normal  to  the  freestream  Ig  is  a  smooth  integral  and 
m  is  unity.  In  general  m  is  a  function  of  the  shock  angle  for  Is 


to  be  a  smooth  function.  IL  involves  only  the  geometry  of  the 
wing.  In  the  formulation  sonic  conditions  exist  when  u  =  1. 
Equation  (126)  can  be  rearranged  to  give 


u  =  [1  i  (1  -  2  m2F)?] 

in 
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and  the  question  arises  as  to  which  root  to  take.  If  the  flow  is 
initially  subsonic  then  u  <  1  and  the  negative  root  must  be 
taken.  As  the  air  flows  over  the  airfoil  the  velocity  increases 
until  u  =  — This  value  is  the  maximum  that  the  negative  root 
in  Equation™ (1 31 )  will  permit  and  occurs  when 


If  the  velocity  is  to  increase  beyond  1/nr  then  the  positive 
root  must  be  taken.  However,  it  is  possible  that  the  velocity 
could  jump  from  a  value  less  than  1/m  by  simply  switching  signs 
in  Equation  (138).  This  would  give  an  expansion  shock  which  is 
not  permitted  since  this  implies  that  entropy  decreases.  Hence 
at  an  "accelerating"  switch,  that  is,  a  switch  in  sign  when  the 
flow  is  accelerating  at  that  point.  Equation  (132)  must  be 
satisfied.  A  compression  jump  is  possible  and  consequently 
Equation  (132)  need  not  be  satisfied  when  the  flow  is 
decelerating.  Finally,  for  the  solution  to  be  real  throughout 
the  domain,  it  is  necessary  that 


(133) 


Equations  (129)  and  (130)  lead  to  the  following  regularity 
conditions  which  eliminate  expansion  shocks: 
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Here,  x  =  xQ,  y  =  y0  is  the  switching  point  for  a  zero  strength 
jump.  If  Equations  (134)  are  satisfied  at  a  decelerating  point 
in  the  flow  then  there  is  no  compression  shock  and  the  flow  is 
shock  free. 


The  transonic  integral  equation  is  given  by  Equation  (126) 
and  if  the  flow  has  a  supersonic  zone  then  the  regularity 
conditions  of  Equation  (134)  must  be  satisfied.  These  regularity 
conditions  in  general  ensure  that  no  discontinuities  are  present 
at  given  switching  lines  xQi,  yQ^  (i  =  1,2).  In  addition  to  the 
u(x,y,  z,  t)  variable  in  equation  (126)  there  are  only  four 
unknowns  for  the  regularity  conditions,  namely  the  value  of  the 
xoi'  Yoi’  Since  six  equations,  (equations  134)  must  be  satisfied 
to  ensure  that  no  discontinuities  occur  it  is  apparent  that  a 
transonic  flow  must  in  general  have  a  discontinuity  because  there 
is  not  a  sufficient  number  of  unknowns  to  eliminate  both 
shockwaves.  This  discontinuity  is  a  compression  shock  if 
expansion  shocks  are  excluded.  However,  it  is  possible  that  the 
necessary  six  regularity  equations  are  satisfied  by  the  four 
unknown  xQ^ ,  y0£ •  In  this  case  the  flow  has  a  supersonic  zone 
but  no  discontinuities.  An  airfoil  with  this  type  of  flow  is  a 
supercritical  shock  free  airfoil. 

Since  u(x,y)  is  a  function  of  the  forcing  term  IL,  Equations 
(134)  can  be  written  as 


1  ▼  „ ( f  +  i  v(iT»  -  0 

Lx  02.  Ol  SX  L  02 ,  02 


I.  (x  . ,y  . ,t)  +  I  (IT,  x  .,  y  . ,t)  =  0 
Ly  02  JQ2  sy  L  02'  J02 


(135  a) 
(135  b) 
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It  (X  <  )  +  I  Ur' 


*~l> 
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(135  c) 
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For  xq  ,y  (i  =  1,2)  these  are  the  six  regularity  equations 
for  thezfou£  switching  points  xq  ,  yQ  and  in  general  the 
equation  set  is  inconsistent.  H&wevei1,  as  noted  above,  it  is 
possible  that  the  flow  is  such  that  Equations  (135)  are  satisfied 
at  both  XQ  '  yo  *  **  interest  to  see  what  happens  if 

such  a  shoci  frei  supercritical  flow  is  perturbed  in  some  way.  A 
perturbation,  either  in  Mach  number  or  geometry,  is  transmitted 
through  the  forcing  term  IL  and  hence  the  effect  of  an 
infintesimal  change  in  1^,  AI^  ,  on  the  solution  is 
investigated.  Since  the  IL  and  Ig  are  functions  that  are 
differentiable  Equations  (135)  can  be  expanded  to  give 

ClLxy(!‘oi,yoi,t)  +  xoi'  yoi’  t):  Mtoi  +  4yoi[ILyy Uoi'  yoi' 


+  I  (I.  ,  x  y  t)]  +  AIT  (x  .  ,  y  t) 
syy  L  oi'  ^oi  Ly  oi  Jox 

"t  Al  Cl*  ,  Al_  ,  X  .  i  >  t)  —  0 
sy  L  L  oi  Ol 


(136  aj 


[I  (x  ,  y  ,  t)  +  I  (I  ,  X  . ,  y  . ,  t)]Ax  +AI  (x  ,  y  ,  t) 
L  O.  O.  S  L  Ol  Ol  O.  L  O.  O. 

XX  1  1  XX  1X11 

+  AI  (I  ,  AI  ,  X  ,  y  ,  t)  + 

S  Li  L  O.O. 

X  11 

AyoiCls  (lL'  Xoi'  yoi'  +  JL  (xoi'  yoi'  t)]  =  0  (136  b) 


.-■.V'.v:. 

V-. 

.  ‘  .1 . • 


AlL(xofyo'.  +  Als(lL'  AIL'  xoi'  yoi'  =  0  (136  c) 

ii 

where  Alg  is  the  change  in  Ig  due  to  the  change  in  IL  of 

AI.  .  In  the  derivation  of  Equation  (136  c) ,  equation  (135)  is 

used.  Equations  (136)  give  six  regularity  conditions  that  must 

be  satisfied  to  retain  the  shock  free  flow,  and  again  only  four 

parameters.  Ax  ,  Ay  (i  =  1,2)  that  can  be  chosen  to  effect 
o.^o. 

these  conditions*  Equation  (136  c)  does  not  contain 
axq  or  AyQ  and  hence  it  is  impossible  to  satisfy  Equation 
(13$  c)  unless  AI ^  is  considered  an  undetermined  function.  Thus 
a  shock  free  flow  cannot  in  general  be  perturbed,  either  by  Mach 
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number  or  geometry,  and  still  remain  shock  free.  However,  if 
AIl  can  be  chosen  such  that  equations  (136  c)  satisfied  for  all 
z  =  constant  lines  covering  the  supersonic  zone  then  Equation 
(136  a,  136  b)  can  be  satisfied  by  a  judicious  choice  of 
Axq  AyQ  ,  thus  it  is  possible  in  certain  circumstances  to 
pertvirb  a  shock-free  flow  and  retain  the  shock-free 
characteristic  if  the  geometry  perturbation  is  suitably  chosen. 
Since  AILis  a  s^n9^e  parameter  function  of  Mach  number  for  angle 
of  attack  (see  Reference  13),  AI^  is  determined  for  all 
(x,y,z,t)  by  a  single  paramter  family.  Thus  it  is  not  possible 
for  a  shock  free  airfoil  to  remain  shock  free  under  a  simple  Mach 
number  or  angle  of  attack  change  since  Equations  (136)  cannot  be 
satisfied  on  all  z  =  constant  lines  where  u  >  1. 


5.2  Summary 


A  theory  for  the  existence  of  shock  free  steady  flow  for 
airfoils  (Ref.  8)  has  been  extended  to  treat  three  dimensional 
low  frequency  unsteady  flow.  The  general  well  known  result,  that 
shock  free  flows  are  singular  points,  unstable  to  small 
perturbations,  is  unaffected  by  the  extension. 

6.0  CONCLUDING  REMARKS 

The  main  objective  of  the  research  reported  here  is  to 
determine  the  relative  magnitude  of  viscosity,  vorticity  and 
entropy  in  a  transonic  calculation.  The  first  conclusion  is 
that,  for  attached  flow,  the  viscous  effect  is  relatively  small; 
the  large  corrections  reported  in  the  literature  being  due  to  the 
presence  of  nonunique  solutions.  Second,  the  vorticity  and 
entropy  effects  are  the  same  order  of  magnitude,  with  the  added 
effect  of  entropy  alone  on  the  shock  strength  and  hence 
location.  The  report  also  considers  the  existence  of  shock  free 
flows  in  three  dimensional,  low  frequency  unsteady  transonic 
flows  and  the  stability  of  such  flows  to  small  perturbations. 
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As  a  result  of  the  investigation  some  insight  into  the 
appearance  of  nonunique  solutions  in  transonic  potential  flow  is 
given.  Also  the  effect  of  entropy  and  vorticity  on  Euler 
equation  solutions  is  studied. 
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APPENDIX 


Further  Comments  On  The  Number  Of  Eigensolutions 


It  is  possible  that  the  only  solution  to  equation  (31)  is 


f(5,n)  =  0 


In  this  case 


Au  u  -  au  —  I  =0 
0  CL 


and 


Au_  = 


■cL 


°  V1 


A'J  =  -  I 


U  -1 
°  J 


/J 


(Al) 


( A2 ) 


(A3) 

(A4) 


If  equation  (A2)  is  substituted  into  equation  (34)  if  follows 
that 


-  i  i  i~Lu)  y 

Au  -  Au  u  =  -  J1  ■ - 

Tf  *  O 


(x-o2+  y2 


( A5 ) 


The  right-hand  side  of  equation  (A5)  can  be  recognized  as 
the  li 
(A5)  becomes 


Au  ,  the  linear  value  of  Au  in  the  flow  field.  Hence  equation 

iw 


AU  ( 1— u  )  =  Au  (A6) 

L 

It  can  be  seen  that  as  H  1  Au  becomes  infinite  which 
will  give  an  expansion  shock.  Consequently,  equation  (Al)  cannot 
be  the  correct  solution  of  equation  (34);  a  nonzero  value  of 
f ( n)  must  exist. 


ITERATION 
SEQUENCE  AND 


CODE 

REMARKS 

GRID 

CONVERGENCE 

CL 

GRUMFOIL 

Viscous 

600/400/500 

0.3  x  10-3 

0.237 

Inviscid 

600/400/500 

0.3  x  10"3 

0.976 

Inviscid 

10/10/800 

0.5  x  10"3 

0.255 

Inviscid 

Viscous  (600- 

400)  Inviscid 

fine  grid,  (500) 

o 

• 

X 

O 

u> 

0.205 

TAIR 

Inviscid 

"C"  grid 

default 

H 

Viscous 

"C"  grid 

default 

0.225 

Inviscid 

"o"  grid 

default 

0.418 

